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1 Introduction 

Polynomial and rational solutions for linear ordinary differential equations can 
be obtained by algorithmic methods. For instance, the maple package DEtools 
provides efficient functions polysols and ratsols to find polynomial and ra- 
tional solutions for a given linear ordinary differential equation with rational 
function coefficients. 

A natural analogue of the notion of linear ordinary differential equation in 
the several variable case is the notion of holonomic system. A holonomic system 
is a system of linear partial differential equations whose characteristic variety is 
middle dimensional. 

Chyzak gave an algorithm to find the rational solutions of holonomic 
systems by using elimination in the ring of differential operators with rational 
function coefficients combined with Abramov's algorithm for rational solutions 
of ordinary differential equations with parameters. To the authors, solving 
holonomic systems is analogous to solving systems of algebraic equations of 
zero-dimensional ideals. Under this analogy, the method of Chyzak corresponds 
to the elimination method for solving systems of algebraic equations. 

The aim of this paper is to give two new algorithms, which are elimination 
free, to find polynomial and rational solutions for a given holonomic system 
associated to a set of linear differential operators in the Weyl algebra 

D = k(2;i, . . . ,a;„,9i, . . . ,9„) 

where k is a sub field of C. 

Polynomial and rational solutions can be obtained, if they exist, by using an 
exhaustive search. For instance, when / = is the singular locus of a holonomic 
system M = D/I, any rational solution has the form g/ ■ If we have upper 
bounds for the degree of the polynomial g and for r, then we can construct all 
rational solutions by solving linear equations satisfied by the coefficients of g. 
Alternatively, if we know the dimension of rational solutions, then we can obtain 
all rational solutions by increasing the degree of g and r. Hence, the problem 
reduces to finding effective bounds for these numbers. 



In sections 2 and 3, we give algorithms for upper bounds on the degree of 
g and on r. The main techniques we use are Grobner deformations in D as 
introduced in the book ijll] and the 6-function for D/I and /. 

In section 4, we give an algorithm to evaluate the dimension of polynomial 
and rational solutions. Our approach is an analog in Z? of a question stud- 
ied by Singer who gave an algorithm to compute IIom/j(M, A'') for left 
R := k(xi)(9i)-modules M and N and studied its relation to factorizations of 
ordinary differential operators. The theory of D-modules translates our problem 
on polynomial and rational solutions to constructions in the ring of differential 
operators D. For example, the k-vector space 

UomoiD/I, k[x]) ~ ®^ T>(D/I)) 

is the space of the polynomial solutions of the left ideal /. Here, is the module 
of the top dimensional differential forms and D is the dualizing functor. See, 
e.g., the book of Bjork ||] on this translation. We evaluate the dimension of 
the right hand side by recent developments of computational algebra such as 
construction of free resolutions in the ring D and restrictions of Z?- modules 
[H' 101 II''' • C)ur method also allows us to evaluate the dimension of solutions 



to a holonomic system inside any holonomic module. For instance, we can find 
the dimension of the delta function solutions to /. 

Throughout the paper, we refer to the book for fundamental facts on 
the algorithmic treatment ol D. Also, the algorithms which appear in the paper 
have been implemented in either kain [|l3| or Macaulay 2 j^. 

We deeply thank Dan Grayson, Anton Leykin, and Mike Stillman, who 
helped to implement D-modules in Macaulay 2, and Frederic Chyzak and Michael 
Singer for discussions on rational solutions. 



2 Polynomial solutions by Grobner deformations 

How can we obtain all polynomial solutions for ordinary differential equations? 
One method is to compute the indicial polynomial at infinity, find an upper 
bound on the degrees of polynomial solutions, and determine the coefficients 
of polynomials. The analogous method works for holonomic systems by using 
Grobner deformations. For £ G D and the weight vector w G R", we denote by 
ini-w,w)i^) the initial term of i with respect to the weight (~w,w) (see, e.g., 
§1.1]). The following proposition follows from the definition of in(_„,^)(£). 

Proposition 2.1 Suppose that f{xi,...,Xn) is a polynomial solution of I = 
D ■ . . . ,£m}- Take w G Z". Then /{f^^xi, . . . ,t"'"a;„) can be expanded as a 
polynomial in t as 

Ui^)tp + 0{tp+'). 

Then we have 

in(-™,t„)(^i) • /«, = 0. 
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The initial ideal in(_^^)(/) is sometimes called the Grobner deformation of 
/ with respect to {—w, w). 

Theorem 2.2 There exist only finitely many Grobner deformations. 

The Newton polytope of a polynomial solution / is defined as the convex 
hull of the exponent vectors of /. For generic w, fw is a monomial cx" and the 
point a is a vertex of the Newton polytope of /. 

Let R — k(xi, . . . , Xn){di, . . . , 9„) and 9i — Xidi. Since a G Z" belongs to 
the zero set of the indicial ideal 

in(-t„, „,)(/) = R ■ \^(-w,w){I) n k[6'i, . . . , Bn], 

we can construct a polytope that contains the Newton polytopes of the polyno- 
mial solutions by taking the convex hull of all the non-negative integral roots of 
all possible indicial ideals. 

It is not necessary to find all Grobner deformations to obtain polynomial 
solutions. Let h{s) be the generator of in(_^,^t„)(J) n k[s], s — Wi^i- The 

polynomial h{s) is called the 6-function of / with respect to {—w,w). The next 
proposition follows from the definition of b{s). 

Proposition 2.3 Let w be a strictly negative weight vector. In other words, we 
assume that Wi < for all i. Consider the b-function b{s) of I with respect to 
{—w,w) and let —ki be the smallest integer root of b{s) = 0. The polynomial 
solutions of I have the form 

Pi>0,P'W<ki 



Algorithm 2.4 (Finding the polynomial solutions by a Grobner deformation) 
Input: a holonomic left ideal /. 
Output: the polynomial solutions of /. 

1. Take a strictly negative weight vector w, compute the Grobner deforma- 
tion in(_„, „,)(/), and compute the smallest non-positive integer root — fci 
of the 6-function with respect to {~w,w). See, e.g., Alg. 5.15] for 
these procedures. 

2. If we do not have such a root, then there is no polynomial solution other 
than 0. 

3. If there is a minimal integer root, then determine the coefficients Cp of (0) 
by solving linear equations for the coefficients. 
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Example 2.5 The following system of differential equations of two variables is 
called the AppcU differential equation Fi(a, b, b', c) iQ: 



Ox{Ox + ey+c-l)- x{0x +0y+ a){ex + b), 

eyiOx + ey + c-i)- y{ex + Oy + a){ey + b'), 

{x - y)dxdy - b'dx + bdy 



where a, 6, 6', c are complex parameters. Let us demonstrate how Algorithm 



2.4 works for the system of parameter values (a, b, b', c) — (2, —3, —2, 5). First, 
we choose a strictly negative weight vector w — (—1,-2) and compute the 
6- function 6(s), s ~ — '29y, which is the generator of the principal ideal 
in(_^^)(/) n Q[— — 2dy]. We can use the V-homogenization or the homoge- 
nized Weyl algebra to get the generator (see, e.g., §1-2]). Second, we need 
to find the integer roots of the fe-function b{s) = 0. In our example, these are 

-7,0,4. 

From Proposition the highest (— u')-degree monomial cx^y'^ in a polynomial 
solution gives rise to an integer solution wip + W2q = — p — 2q of the ^-function. 
Hence, the polynomial solutions are of the form 

p,q>a,p+2q<7 

Finally, we determine the coefhcients Cpq by applying the differential operators 
to / and putting the results to 0. In our example, we have only one polynomial 
solution 

/ 1 2 1 4 , 3 , 3 2 24 3,2 

( y + -y )x^ + ( — y y + -)x 

^ 2r V Zh' 44^ 35^ W 

/ 12 2 6 6, 1,4 
+ {-^^y'^-^y-^)x^-y---y^X. 



3 Rational solutions by Grobner deformations 

The singular locus of a D-ideal / is defined to be the projection of the char- 
acteristic variety of / minus the zero section from the cotangent bundle to the 
coordinate base space. In other words, it is the zero set 

Sing(/) =y((in(o,,)(/) : (6, ■ • • , Cn)"") H k[xi, . . . , a;„]) . 

Any rational solution to / has its poles contained inside the singular locus. 
Thus if /(x) defines the codimension 1 component of Sing(/), we may limit our 
search for rational solutions to k[x][j]. 

We will present a method to obtain an upper bound of the order of the poles 
along f — for each rational solution. For this purpose we use the notion of the 



4 



6-function for / and a section u of a holonomic system, which was introduced by 
Kashiwara : Let V be the sheaf of algebraic differential operators on X = C" . 
For a holonomic P- module M = V/VI and a polynomial /, consider the tensor 
product 

M^o[r\s]r<^o^M. (2) 

This N has a structure of a left P-module via the Leibnitz rule. Let u be a 
section of Ai. Then the ^-function for / and u (or for f^u) at p S C" is the 
minimum degree monic polynomial ^ b{s) G C[s] such that 

(3) 

holds in A/" at p (i.e., as a germ of A/" at p). This fe- function depends on the point 
p. As a function of p, there is a stratification of C" for which the 6- function does 
not change on each strata (see e.g. Q for an algorithmic proof of this fact). In 
the definitions (^) and (||) for 6-function, if we replace O by the polynomial ring 
k[x], T) by the Weyl algebra Z), and Al by a holonomic Z3-module M — D/I, 
then we obtain the global 6-function for / and u. It is the least common multiple 
of 5- functions at every point. 

Theorem 3.1 Let u be the residue class of 1 in V/VI , and let 6(s) be the b- 
function for f and u at a point p G C" where f{p) = 0. Assume that I admits 
an analytic solution of the form gf^ around p, where r d C, g is a holomorphic 
function on a neighborhood of p, and g(p) ^ 0. Then s + r + 1 divides b{s). 

Proof : Let I?^" and be respectively the sheaf of analytic differential opera- 
tors and the sheaf of holomorphic functions on C". We may define the analytic 
6-function by replacing O by C'^, V by V^'^ and 7W by a I?™-module in 
the definitions ^ and (^). Since the 6- function is an analytic invariant and the 
analytic and the algebraic 5- functions coincide (see e.g. H, §8]), we may work 
in the analytic category. We do this to consider solutions gf^ where g is holo- 
morphic at p. If we only wish to consider solutions gf^ where g is a polynomial, 
then we may work in the algebraic category. 

In general, given a map of left ^'^"-modules : A^f" A^f" and a section 
u of Aif^, the 6- function for f^u at a point p is divisible by the 6- function for 
f^4'{u) at p. We apply this basic fact to the following map Lp. Let J™ be the 
annihilating ideal of gf in P'^". Since J''" D J'*" := and g{p) ^ 0, we 

have a left I>'^"-homomorphism 

which sends u to gf^. This map extends to a left I?''"[s]-homomorphism 

1(E) if : 0^"[/-i,s]/'* ®oa„ P'^"//^" 

which sends f^ E)u to gf^^ . By the definition of 6(s), there exists a germ P{s) 
of V[s\ at p such that 

{P{s)f-b{s)){r®u) = Q. 
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Since 1 (8) (/? is a left X>™-homoniorphism, applying it to the above equation gives 
the equation {P{s)f — b{s)){gf''~^^) = 0, or in other words, 



Thus, we see that the Bernstein-Sato polynomial bf(s) of / at p divides b{s — r). 
Note that s + 1 divides 6/(s) since f{p) = (cf. Q). In conclusion, we have 
proved that s + 1 divides b{s — r). This completes the proof. □ 

By virtue of the above theorem, we can obtain upper bounds by computing 
the fe-function for f^u at a smooth point of each irreducible component of the 
singular locus of /. From now on, let us also take / e k[x] to be a square-free 
polynomial defining the codimension one component of the singular locus, and 
let f — fi ■ ■ ■ fm be its irreducible decomposition in k[x]. 

Theorem 3.2 Let bi{s) be the b-function for ffu at a generic point of fi = 0. 
Denote by Vi the maximum integer root ofbi{s) = 0. Then any rational solution 
(if any) to I can be written in the form gfi^^~^ ■ ■ ■ fm^"^^^ with a polynomial 
g G C[x]. If some bi{s) has no integral root, then there exist no rational solutions 
to I other than zero. 

Proof: An arbitrary rational solution to M is written in the form gfi'^^ ■ ■ ■ fm""" 
with integers z/i , . . . , Vm and g € C [x] . Since the space of the rational solutions 
with coefficients in C is spanned by those with coefficients in k, we may assume 
g £ k[a;], and / and g are relatively prime in k[x]. Let p be a generic point of 
fi — 0. We may assume that fi is smooth at p, g{p) ^ 0, and fj{p) ^ for 
j ^ i. It follows from Theorem that bi{vi — \) = This implies Vi < ri + 1. 
□ 

Since fe-functions divide the global 6-function, an upper bound can also be 
obtained from the global fe-function. 

Corollary 3.3 Let bi{s) be the global b-function for ffU, and denote by ri the 
maximum integer root ofbi{s) — 0. Then any rational solution (if any) to I can 
be written in the form gfi^^^^ ■ ■ ■ fm^""~^ with a polynomial g G C[x]. 

We mention the corollary since the algorithm to compute global 6-functions 
is simpler than the algorithm to compute 6-functions. However, the ^-function 
offers finer information. For instance, the well-known example f = + y'^ + 
z^+uo^ has Bernstein-Sato polynomial (s-|-l)(s-|-2) coming from the functional 
equation, \{dl + dl + dl + dl) ■ = (s + l)(s + 2)/". Now consider the 
module M = D ■ f^^ and let u be the section of /^^. The global 6-function for 



f^u is s{s + 1), and hence Corollary 3.3 implies that rational solutions of M all 



have the form gf^^ or g/*^, where 5 is a polynomial not divisible by /. On the 
other hand, the Bernstein-Sato polynomial of / at any nonsingular point p of 
/ = (i.e. except for the origin) is s -I- 1. It follows that the 6-function for 



equals s at the generic point of / = and hence Theorem 3.2 implies that all 
rational solutions actually have the form gf^^. 
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An algorithm to compute the 6- function and the global &- function for 
was first given in ||] based upon tensor product computation, which is slow and 
memory intensive. Shortly thereafter, Walther introduced in a more efhcient 
method to compute the global 5-function for fu. Both methods give the global 
5-function exactly, under the condition that / is /-saturated. Otherwise, we get a 
multiple of the global &- function. Similarly, the method of Q gives the &- function 
exactly if / is /-saturated and additionally a certain primary decomposition in 
C[x] is known. If primary decomposition is only available in k[x], we again get 
a multiple of the 6-function. 

Let us now describe an algorithm to compute the 6-function for /"u at 
a generic point of / = by combining the method of [p^ and the primary 
decomposition as was used in 

Algorithm 3.4 (Computing an upper bound of the 6-function at a generic 
point) 

Input: a finite set Gq of generators of a holonomic Z3-ideal / and an irreducible 
polynomial / G k[x]. 

Output: 6'(s) e k[s], which is a multiple of the 6-function 6(s) for f^u at a 
generic point of / = 0, where u is the residue class of 1 in D/I. 

1. Introducing a new variable t, put 'di = di + [df /dxijdt- Let / be the 
left ideal of Un+i, the Weyl algebra on the variables xi, . . . , a;„, i, that is 
generated by 

{P(x, . . . , z?„) I P(x, 9i, . . . , a„) e Go} U {i - /(x)}. 

2. Let Gi be a finite set of generators of the left ideal in(_i^o,...,0;i.o,...,o)(^) 
of Dn+i- Here, —1 is the weight for t and 1 is the weight for dt- 

3. Rewrite each element P of Gi in the form 

P = d^P'{tdt,^,di,...,dn) or p = t^p'{tdt,x,di,...,dn) 
with a non- negative integer fj,, and define ipiP) by, 

V'(P) := i^af P' ^tdf- (tdt - M + l)P'itdt, X, 9i, . . . , 9„) 
or ^{P) := d^t^'P' = {tdt + 1) • • • {tdt + ^l)P'{t^t,^. . . . , 

Put 

G2 :={7^(P)(-s-l,x,9i,...,a„) l^'GGi}. 

4. Compute the elimination ideal J := k[s,x] n D[s\G2- (The global 6- 
function can be obtained at this stage by computing the monic generator 
of the ideal J n k[s].) 

5. Compute a primary decomposition of J in k[s,x] as 

j = Qin---nQ^. 
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6. For each i = 1, . . . , v, compute Qix '■= Qi H k;[x], which is a primary ideal 
of k[x]. 

7. Let b'{s) be the monic generator of the ideal 

f]{Q^ n k[s] I y/Q~ C k[x]/} 

of k[s]. (Note that k[x]/ implies that \/Qix equals k[x]/ or {0}.) 

Theorem 3.5 In the above algorithm, the polynomial b'{s) is precisely the b- 
function for f^u at a generic point of f = if I is f -saturated (i.e., I : f°° — I) 
and each C[s,x]Qi remains primary in C[s,x]. Otherwise, the polynomial b'{s) 
is a multiple of the b-function for fu at a generic point of f = 0. 

Proof: Using essentially the same method as the proof of Lemma 4.1 in JlGf , 
we can prove that / is precisely the annihilator ideal for 6{t — /(x)) tt in 

M := - /(x))) ®c[x] D/I, 

where S{t—f{x)) denotes the residue class of (<— /(x))^^ in k[x, (f— /(x))^-^]/k[x]. 
Let bt{s) be the indicial polynomial for S{t — /(x)) (g) u along < = at a point 
{0,p) with fip) = 0. Then by Theorem 6.14 of @, the 6-function b{s) for fu 
at p divides out, and if / is /-saturated, coincides with bt{—s — 1). 

It follows from the definition that bt{—s — 1) is a generator of the ideal 
Op[s]J n C[s] of C[s], where Op denotes the stalk of O at p. If C[s,x](5i are 
primary in C[s,x], fe'(s) generates the above ideal in view of Theorem 4.7 of 
[|| (cf. also Lemma 4.4 of Q). In general, although Qi is primary in k[s,x], 
the extension C[s,x](5i is no longer primary in C[s,x] and admits a primary 
decomposition 

c[s,x]Q, = Q,in---nQ,^,. 

In this case, b'{s) is the least common multiple of the generators of the ideals 
Op[s]Qij n C[s] for j — 1, . . . , ^i, while 6t(— s — 1) is the generator of Op[s]Qij (1 
C[s] for some j (such that p belongs to the zero set of Qij C[x] which is also 
the zero set of a factor of /). This completes the proof. □ 

Remark 3.6 In the notation in the above proof, the linear factors of b'{s) and 
those of bt{—s — 1) coincide. In particular the set of integer roots of b'{s) = is 
the same as that of bt{—s — 1) = 0. In fact, this follows from the fact that the 
linear factors of fo'(s) in k[s] are invariant under the action of the Galois group 
of k over k. 

Remark 3.7^is /-saturated if and only if the —1-th cohomology group of the 
restriction of M in the proof of Proposition |3.5| to t = vanishes (see Theorem 
6.4 and Proposition 6.13 of Q), which is computable by Algorithm 5.10 of Q 
or by Algorithm 5.4 of 
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Remark 3.8 The /-saturation of /, which is the ideal D[j] ■ I H D, may be 
computed by using the locahzation algorithm of |10(| (if I is specializable along 
/) or by using the less efhcient algorithm of jlj] (if / is general). By replac- 
ing / with its /-saturation, we may then compute the local ^-function exactly. 
However, since saturation is often an expensive algorithm, we avoid making this 
replacement in practice. 



Once we have determined the integers ri, . . . , of Theorem 3.2, we can use 



Grobner deformations to obtain the rational solutions. Put ki = r^ + l. Then by 



virtue of Theorem 3.2, we have only to determine rational solutions of the form 
gfi'^^ ■ ■ ■ fm^"^ f'^^ some polynomial g. This amounts to computing polynomial 
solutions of some twisted ideal /(fci,. of /. Namely, consider how fdi acts 
on an clement gfi^^ ■ ■ ■ fm'^'" ■ 



In other words, fdi acts on the numerator g as the differential operator 

' fj dx 



Thus, if we multiply a set of generators {171, ... ,3™} of / by sufficiently high 
powers /™' such that /™'5i G k(xi, . . . , a;„, fdi, . . . , /9n), and if /' is the ideal 
obtained from {/™^(?i, . . . , f™"'gr} by substituting Li for fdi, then the rational 
solutions of / have the polynomial solutions of /' as numerators. Thus, it only 
remains to compute the polynomial solutions of /'. However, since /' might not 
be holonomic, we cannot apply Algorithm |2.4| just yet. 

Hence let us define I{kx,...,km) •= k(x)(9) -I' C\D, which is the Weyl closure of 
/' and whose polynomial solutions are the same as /. The advantage of /(/^^ /;^) 
is that it is indeed holonomic, which follows from a theorem of Kashiwara. 
Namely, note that since / is of finite rank, /' remains of finite rank because an 
element in / of the form {gi{x)df' + lower order elements) will be sent to an 
element in /' of the form [f^^gi{'x)df^ + lower order elements). Now let ft,(x) be 
any polynomial vanishing on the singular locus of /'. Then the non- holonomic 
locus of /' is contained inside the zero set of /i(x) regarded as a function on 
the cotangent bundle, and a theorem due to Kashiwara states that the ideal 
D[h~^] ■ I' O D is holonomic. Furthermore, an argument in shows that the 
Weyl closure of /' also equals D[h~^] ■ I' n D. Summing up, we arrive at the 
following algorithm. 

Algorithm 3.9 (Computing the rational solutions of a holonomic ideal) 

Input: generators of a holonomic D-ideal /. 

Output: A basis of the rational solutions h e k(x) of / • ft, 0. 

1. Compute a polynomial / defining the codimension 1 component of Sing(J). 
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2. Compute the irreducible decomposition ,f = fi ■ ■ ■ fm in k[x]. 



3. For each i — 1, . . . ,m, compute the output b'{s) of Algorithm 3.4 with / 
and fi as input. Let be the maximum integer root of b'{s) = and put 
ki = ri + 1. If b'{s) has no integral root for some i, then there exists no 
rational solution other than zero. 

4. Compute the twisted ideal I(ki....,km) follows. First, form the ideal 
/' described in the paragraphs preceding the algorithm. Second, compute 
any polynomial /i(x) vanishing on the singular locus of /'. Third, compute 
the localization {D/ I')[h^^] using the algorithm in [|loj. Then the ideal 
I(ku-M is the kernel of D ^ D/I' iD/I')[h-^]. 

5. Compute a basi s {g i, . . . , g^} of the polynomial solutions of I{ki,...,k^) 



using Algorithm 2.4 



6. Output: {gi/i J^^"^ gkfi a basis of the rational 

solutions of /. 



Example 3.10 Let / be the left ideal generated by 

Li = e^{e^ + Oy) - xi9, + ey + 3)(0, - 1) 

and 

L2 = ey{e, + By) - y{e, + 6y + i){ey + 1). 

The Appell function Fi (3, — 1, 1, 1; a;, is a solution of this system. The singular 
locus of / is xy{x — y){l — x){l — y) = 0. We can compute the local indicial 
polynomial of u, the modulo class of 1 in D2/I, along a; = directly by the 
algorithm of |||, Section 4]: It is s{s — 1) on {{0,y) | y ^ 0}, and s{s — 1)^ at 
(0,0). In the same way, the indicial polynomial of u along ?; = is s(s + 1) on 
{(a;,0) \ x^0}, ands(s + l)(s-l) at (0,0). 

Now let us compute the 6- function for (l—y)'^u. The local indicial polynomial 
oi d{t + y — 1) IS) u along t = is s(s + 3) at any point of t = 0. Hence the 
6-function for (1 — yYu divides (s + l)(s — 2). In the same way, the local indicial 
polynomial of 6{t + x — 1) (g) u along t = is s{s + 1) at any point of t = 0. 
Finally, the indicial polynomial of S{t — x + y) u is s(s — 1) on {{x,x) \ x ^ 0}, 
ands(s-l)(s-2) at (0,0). 

Therefore, we conclude that any rational solution to /, if it exists, can be 
written in the form g(x, y)y~^{l — a;)~^(l — y)~'^ with a polynomial g. Now we 
may compute the twisted ideal /(o, 1,0, 1,3) j where fi = x, f2 = y, fs = x — y, fi = 
X— 1, /5 ~ 1, and / is the product. Multiplying by we get the expressions, 

f^L, = (x2-x3)(/a,)2+x((l-3x)/-(l-x)yg-(l-x)xfi)(/9,) + 

x{l - x)y{fdy){fd,) + xyfifdy) + 3x/2 
f^L2 = {y'-y'Kfdyr + y{{l-5y)f-{l-y)x^-{l-y)y%){fdy)+ 

y(l - y)x{fd,)ifdy) ~ yxfifd.,) - 3yf^ 
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and we set Ti and T2 to be the operators obtained from the substitution (Q). 
We remark that the ideal generated by Ti and T2 is neither holonomic nor 
speciahzable with respect to the weight vector (1, 1, —1, —1), so it is difficult to 
apply Grobner deformations to it just yet. 

The twisted ideal /(o, 1,0. 1,3) is the Weyl closure of the ideal generated by Ti 
and T2. At the moment, this Weyl closure is computationally too intensive to 
compute. However, we are able to compute a partial closure by noting that Ti 
is divisible hy g = /2/I/4/5 and T2 is divisible by h = /i /|/|/4/5- Now the 
ideal J generated by 

iTi = {-x^y + + 2x'^y ~ 2x'^ - xy + x)dl + 

(-x2y2 + x^y + 2xy^ - 2xy - + y)d^dy + 
{3x^y - 6xy + Zy)d, + {2xy^ - 2xy - 2y^ + 2y)dy+ 
(-4xy -2x + Ay + 2) 
= i-xy'^ + 2xy^ + y^~xy^-2y^ + y^)d^ + 

{—x^y^ + 2x'^y'^ + xy^ — x'^y — 2xy'^ + xy)dxdy+ 
{ix'^y'^ — Ax'^y — 3xy^ + + Axy — x)dx + 
(Axy^ - &xy'^ - Sy^ + 2xy + - y)dy+ 
{-6xy^ + 8xy + 3y^ - 2x - 4y + 1), 



is indeed holonomic, hence we may apply Algorithm 2.4. We find that the b- 
function with respect to the weight w = (—1, —1) is (s + 5)(s + 2)^(s + l)s(s — 
2)(s — 3)^, which implies that a polynomial solution to J must have degree less 
than or equal to 5. We find that J has 2 polynomial solutions, so that / has 2 
rational solutions, {xy"^ — 3xy + 3a: — l)/{y — 1)^ and {x — y)/y{y — l)"^. 

4 Solutions by duality 

For holonomic M and N , it is well known |0] that 

Extl,(Af, N) ®^ (D(M) iV)), (5) 

where {D/{xi^ . . . ,x„} • D), and D(M) is the holonomic dual, 
D(M) :=Homk[,](r!,ExtS(Af,i?)). 

The spaces Ext^ (Af , N) arc finite-dimensional k- vector spaces and corre- 
spond to the solutions of M in N when i — 0. For example, if TV = k[x], then 
we obtain the polynomial solutions of M, whereas if iV = D/D ■ {xi, . . . 
then we obtain the delta function solutions of M with support at the origin. 

In this section, we explain how (||) can be used to compute the dimensions of 
Ext^ (Af, A^) . We first discuss how to compute the holonomic dual, next discuss 
the special cases A^ = C[x] and A^ = C[x][j], and last discuss the general case of 
holonomic N . A method to extend these algorithms to compute an explicit basis 
of Hom£)(Af, A^) and Ext^(Af, A^) is the subject of the forthcoming paper |p^ . 

Notation: Let us explain the notation we will use to write maps of left or 
right Z?-modules. As usual, maps between finitely generated modules will be 
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represented by matrices, but some care has to be given to the order in which 
elements are multiphed due to the noncommutativity of D. 

Given an r x s matrix A = [aij] with entries in £), we get a map of free left 
D-modules, 

D-^-^D" \gi,...,gr\^[gi,...,gr\- A, 
where and are regarded as modules of row vectors, and the map is matrix 

■A 

multiplication. Under this convention, the composition of maps — > and 

S AB 

map D"^ — > £)* where AB is usual matrix multiplication. In 

general, suppose M and N arc left _D-modules with presentations /Mq and 
/Nq. Then the matrix A induces a left D-module map between M and N, 

denoted (D^'/Mq) {D^ /No), precisely when g ■ A e Nq for all row vectors 
g e Mq. Conversely, any map of left i'-modules between M and N can be 
represented by some matrix A in the manner above. 

Now let us discuss maps of right £)-modules. The rxs matrix A also defines 
a map of right iJ-modules in the opposite direction, 

{dY ^ (dT [hu hsf ^ A ■ [/ii, . . . , hsf, 

where the superscript-T means to regard the free modules (f*)"^ and {D^)'^ as 
consisting of column vectors. This map is equivalent to the map obtained by 
applying Homi)(— ,£)) to D"^ — > D*, thus {D^)'^ may be regarded as the dual 
module Hom£)(£'*, -D). We will suppress the superscript-T when the context is 
clear. As before, the matrix A induces a right iP-module map between right D- 
modules N' = (r>*)^/iV^ and M' = {D'-)'^/M(, when A • 5 e for all column 

vectors g€ N^^. We denote the map by {D'Y (£)'')^/M^. 

Left-right correspondence and fi: As is well-known, a standard use for 

J7 is to establish a correspondence between the categories of left and right D- 
modules. The correspondence can be expressed through the adjoint operator t, 
which is the algebra involution 

t:D — x°'d'^ ^ {-dfx". 

Namely, given a left D-module M ~ D'^/Mq, the corresponding right £)-module 
is f2®k[x] M ~ D^/t{Mo). Conversely, given a right £)-module A'' c± D^/Nq, the 
corresponding left D-module is Homi(.[x](i^, A^) — D''/t{No). Similarly, given a 
homomorphism of left i3-modules (p '■ {D^ / Mq)^(D^ /Nq) defined by left multi- 
plication by the rxs matrix A = [uij], the corresponding homomorphism of right 
D-modules t((/.) : (I?'^/t(Mo))^(£'Vt(Ao)) is defined by right multiplication 
by the s x r matrix t{A) := [T{aij)]^ . 

Let us explain details of the above correspondence for the non-specialist. 
Given a left _D-module M, there is a corresponding right Z)-module fl (S)k[x] M 
where the structure is given by extending the actions, 

{w (g) m)f = wf (g) m (w ® m)$, = w^^m — w ^ 
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for / € k[x] and ^ G Dcr(k[x]). Given a presentation D'^/Mq for M with 
generators denoted {ei}[^i, then in Q (8)k[x] -^'-^ wc have 

(1 ei)x"d'^ = (1 O a;"ei)a'^ = 1 O (-a)^a:"ei = 1 O r(a;"a'^)ei. 

It follows that fl (S)k[x] is generated by {1 (g) eiW^i and gets the presentation 
^yr(Mo). 

Conversely, given a right f-module A'', there is a corresponding left D- 
module Homi5;[x](f^, where the structure is given by extending the action, 

{fip){w) = ip{w)f {£,^){w) = ^{wO - 

for ifi G Homk[x](f^, -^), w £ f G k[x], and ^ G Der(k[x]). A morphism 
(f G Honik[x](0, iV) can be identified with its image <p(l) G N. Since 

= ip{i)T{x-d0), 

the morphism x^d^cp gets identified with (p{l)T{x"d^). In particular, given 
a presentation D'^/Nq of N, then IIomk[x](ri, A^) is generated as a left D- 
module by the morphisms {y'ijl^i such that = ej. By the computation 

above, a relation ^ - agi = in N corresponds to a relation J2i '^{9i)fi in 
Homk[x](^^, A^) because {J2i'ri9i)'Pi){^) = = Yl,i^i9i- It follows 

that IIomk[x](f^, AT) is generated by {</?i}|=i and gets the presentation 

Homk[x](J^,A^)^i5:/T(Aro). (6) 

4.1 Holonomic dual 

Let us discuss how D(M) can be computed. 

Algorithm 4.1 [Computing the holonomic dual] 

Input: D^^ jD ■ {^i, . . . ,5^1}) a presentation of a holonomic left D-module M. 
Output: The holonomic dual D(M). 

1. Compute the first n + 1 steps of any free resolution of M. Let the n-th 
part of the resolution be Dp ^ D" ^ D"-. 

2. Dualize and apply the adjoint operator (recall if P = [pij], then r(P) = 
[T{pij)f) to get DP D« -l^^ D\ 

3. Return ker(-T(P))/Im(-T(Q)). 

Proof: Let the first n + 1 steps of a free resolution of M be denoted, 

F' : D'^"-! ^ > D'^o ^ 0. 
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Applying HomD(D,— ) yields a complex of right D-moduIes, 

RomD{D,F') : {D'-^+^f ^ {D'--f ^ {D'-"-'f < ^ {D'^'f ^ 0, 

and by definition, 

Since D(Af) = Homi(.[x](f^, Ext^(M, D)), it only remains to determine the ef- 
fect of applying Homi(.[x](^^, — )• Using the equation if {Li, . . . , Lk} are gen- 
erators of if = ker(D''"+i D''"), and J2i Ugi € / = Im(i:)''" 4^ Z?''--!) is a 
relation, then the corresponding relation T{gi)ipi in Homi(.[x](f^, Ext2)(M, D)) 
can be realized as the relation J2i'''i-^i9i) — '''{9i)'^{^i) ^ '^{^)- It follows that 

D(M) ~ ^ 



which is the output of step 3. □ 



Example 4.2 The Appell differential equation Fi{2, —3, —2, 5) of Example 
has the resolution O-^ZJi L>2 -Qo jji _^ q ^j^^^.^ 



{6^ - 3)dy - {9y ~ 2)5, 

(2/2 - -I- d^) - 2{y + x)d,, + Aydy + 29, 



(y^ - v){dxdy + dl) - 2xdx + 6ydy + (9, - 9dy 

-{9., ~ 5)dy + i0y ~ l)d^ 

The holonomic dual D(Fi(2, — 3, — 2, 5)) is the cokernel of t{Qi) and is the 
Appell differential equation 1, 4, 2, —3). 



4.2 Polynomial and rational solutions by duality 

When N ~ k[x], the isomorphism (H) specializes to 

Ext^(A/, k[x]) ~ H'^-'ifl (g)^ D(M)). 



(7) 



The right hand side is equivalently the {n — i)-ih integration of D(M) to the 
origin. An algorithm to compute integration is given in Using it, we can 
evaluate the dimensions of Ext^(Af, k[x]) and in particular HomD(A/, k[x]). 

Algorithm 4.3 [Evaluating dimensions of polynomial solution spaces] 
Input: a holonomic left -D-module M. 
Output: dimensions of Ext^(Af, k[x]). 
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1. Compute the dual D(M) using Algorithm 4.1 



2. Compute the integrations of D(Af) to the origin using the algorithm in Q . 
They are finite dimensional vector spaces. 

3. Return the dimensions. 

The dimensions of rational solution spaces can be evaluated in a similar way. 
When N = k[x][j], the isomorphism (||) specializes to 

Extl,(M,k[x][l//]) ^ H--\n^^r, mm/f])- (8) 

The right hand side is now equivalently the (n — i)-th integration of D(M)[j] to 
the origin. An algorithm to compute localization is given in [|lo). Using it and 
the integration algorithm, we can evaluate the dimensions of Ext^(Af, k[x][j]) 
and HomD(A/,k[x][i]). To get the dimension of all rational solutions, take / 
to be any polynomial vanishing on the singular locus. 

We summarize how to compute the integration of a module N to the ori- 
gin according to in a slightly more general way. The generalization some- 
times gives a more efficient strategy than Q. We need to recall some defini- 
tions. To any strictly positive w g Z"q, we get an integration filtration Fj^, 
of D defined by F^{D) = Spanj,{a;"9''|w ■ a — w ■ (3 < i}. More generally, 
for fh G Z"^, we also get a shifted filtration F^[to] of the free module de- 
fined by F:^[rh]{D^) — SpeLni^{x"d^ej\w ■ a — w ■ P — mj < i}. We wiU often 
write D'"[m] for the free module equipped with the shifted filtration F^[m] 
when the context is clear. The filtrations i^u,[m] induce filtrations on subquo- 
tients of in the natural way. Now we may say the steps of the integration 
algorithm. First, compute a (w, — w)-strict free resolution G* of TV of length 
n + I. This is a resolution of N by free modules D^^[mj] with the property 
that the differentials preserve the filtration and moreover induce a resolution 
on the associated graded level. Second, compute the integration 6-function of 
TV with respect to {w, —w), and find its minimal and maximal integral roots fcg 
and ki. The integration 6-function is the monic polynomial 6(s) of least degree 
satisfying h{Y^^WidiXi) ■ F'^{N) C F^^{N). Third, compute the cohomology of 
the complex F-koi^ G')/F^ki-i{^®D G*), which is a complex of finite- 
dimensional vector spaces. The dimensions of the cohomology groups are equal 
to the dimensions of the integration modules of N . 

Example 4.4 Let us evaluate the dimension of polynomial solutions to the 



Appell differential equation M = -^1(2, —3, —2, 5) of Ex ample 2.5. Choose the 
weight vector w = (1, 2). The resolution of Example [4.2| , after dualizing, apply- 
ing the adjoint operator, and shifting, 

0^£)i[0] -^^^ D^[-l, 1] '^^^ D'[0] ^ 0, 

preserves filtrations but does not induce a resolution on the associated graded 
level. On the other hand, if we adjust the resolution to 

G' : 0^D^[1] ^ D^[0, 1] ^ D^[0] ^ 0, 
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where 



Pn 



Pi 



-{e, + ^)dy + {ey + 2)d:, 

{x^ - x){dl + d^dy) + Axd-^ + 2(3.x + 2y)dy + 49^ - 55^, 
(a;2 - x){dl + a.^j,) + 2xd^ + 4(a; + j;)^^ + 59, - Ady - 

{9^ + A)dy - [Oy + 2)9, 



then we do obtain a (w, — ■u;)-strict resolution of D(M) = 1, 4, 2, —3). 

The integration 6-function with respect to {w, —w) is (s + 5)(s — 2)(s — 5), 
hence the integration complex for D(i\f) is quasi-isoniorphic to the truncated 
complex F^(ri(g)£i G*) / {Vl ® n G*), which is a complex of finite-dimensional 
vector spaces with dimensions, 



0^ Q 
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^ ^ 



For instance, F^(r2[0]) consists of the 12 monomials, 

{1, X, y, x^, xy, y^, x^, x'^y, xy"^ , x'^, xy^, x^}, 

and so on. Note that r(Pi) is a {w, — w)-Gr6bner basis of Pi (2, —3, —2, 5) and 
hence for this case, the duality method essentially coincides with the Grobner 
deformation method of Section ^ at the level of Homo(M, k[x]). The above 
computations were made in Macaulay 2, where we get the output. 



11 : PolyExt(M) 
ol = HashTable{ 



=> 

1 => 

2 => 



Here, the output i => qqJ means that dimExt^(M, k[x]) ~ j. 



Example 4.5 Let us now evaluate the dimension of rational solutions to M = 
Pi (2, —3, —2, 5). The singular locus is xy{x — y){x — l){y — 1). We will search 
for solutions in k[a;,y][i] first. From Example [4.2| , D(M) has the presentation 
D/t[Qi). Let u be the section corresponding to the residue class of 1 in this 
presentation. Then the localization D(M)[i] is generated hy ^ and gets 
the presentation D / J , where 



32 



20, + 12) - (0, +ey+ A)dy 



\ {ejy + 20, + lOy + 14) - + \{))xdy 

The natural localization map can be written as ip : D/t{Qi) — > D/ J, where 
</5(l) = x"^. Choose the integration weight vector w — (1, 2). Then D(M)[i] has 
a (w, — ■u;)-strict resolution 

G" : -> ^'[0,-1] ''"^'^ PMO] - 
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where 

ui = —x'^dxdy + xydxdy + 2xdx — llxdy + 7ydy + 14: 

U2 = x^dx + x^dxdy — x'^d'l — x'^dxdy + l&x'^dx + llx^dy+ 

4:xydy — 9xdx — llxdy + 52x — 7 
vi = x]^dl + x^^dxdy - x'^dl - x'^dxdy + Wx'^dx + 12x'^dy+ 

4xydy — 8xdx — llxdy + 52x — 6 
V2 = x'^dxdy — xydxdy — 2xdx + llxdy — Qydy — 12 

The integration 6-function is (s+12)(s+5)(s+2), hence we want the cohomology 
of the complex (g)£) G')/F^{Q. (g)D G'), which has the shape, 

^ Q41 _^ q86 _^ q45 ^ 0^ 

By evaluating the dimensions of the cohomology groups in Macaulay 2, we find 

12 : RatlExt(M,x) 

02 = HashTable{ => QQ^ } 

1 => qq^ 

2 => qq^ 

Since wc already computed a polynomial solution, this means there is one ra- 
tional solution with pole along x. Similarly, we get the exact same dimensions 
for Ext^(M, k[a;, which means that there is also one rational solution 

with pole along y. The rank of the system is 3, therefore we have found all the 
solutions. We could also compute, 

13 : RatlExt(M,f) 

03 = HasliTable{ => qq' } 

1 => qq^ 

2 => qq= 

where / is any of the polynomials x — y, x — 1, or y — 1. As expected, there 
are no rational solutions with poles along x — y, x — 1, or y — 1, but in all cases 
there are new Ext^ and Ext^. We have not computed Ext with respect to any 
products of poles since it is computationally too intensive for now. 

Once we have evaluated the dimension of the solution spaces, we can compute 
the solutions by a brute force method. 

1. For a given holonomic system M, compute its singidar locus. Let / be a 
polynomial such that / = contains the singular locus. 

2. Evaluate the dimension d of the rational solutions by the homological 
duality method. 

3. Try to find rational solutions of the form r{x)/ f'^ , degree(r) = p. Increase 
p + k until we find d linearly independent solutions. 
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4.3 Holonomic solutions by duality 

The isomorphism (||) can also be expressed (see e.g. as 

Extl3(A/, N) ~ Tor^_,(ExtS(Af, D),N). 

To compute the right-hand-side, it is also well known that for all Af, and in 
particular M' = Ext2)(Af, D), 

Tor^_,(Af', N) ~ H\K'{{M' ®k[x] n-^)®N- {x, - y„ 8, + (5.}r=i)), 

where if* denotes the Koszul complex and ® denotes the external tensor prod- 
uct into the category of D2n = k(a;i, . . . , ?/„, 9i, (5i, . . . , i9„, i5„)-modules. 
Combining these isomorphisms leads to 

Extl,(Af, N) ~ ff'(if((D(M)SiV;{x, -2/„a, +<5jr=i)) 

By an automorphism of 13, we can transform {xi — j/i, + into {xi, yi\^^i, 

for which the Koszul complex computes the derived restriction to the origin. 

Algorithm 4.6 [Evaluating dimensions of holonomic solution spaces] 
Input: holonomic left D-modules M and N 
Output: dimensions of Ext^ (Af, A^) . 

1. Compute the dual D(Af) using Algorithm 4.1 

2. Form the Z32n-niodule D(Af) (g)^ N and apply the change of coordinates 
'7 : — > whcrc 77 maps, 

^ \xi - 5i, di i-> + 

3. Compute the restrictions of 77(D(Af) (g)fe A) to the origin using the algo- 
rithm in 1^. They are finite dimensional vector spaces. 

4. Return the dimensions. 



Exampl e 4 .7 Let M — i^i(2, — 3, — 2, 5) be the Appell differential equation of 
example 2^, and let A — kfx, y][i]/k[a;, y]. It has presentation DjD ■ {x^dy\, 



where the generator 1 corresponds to ^ . Using the above algorithm, we compute 

14 : DExt(M,N) 

04 = HashTable{ => QQ' } 

1 => qq' 

2 => QQ= 

Similarly, let A — \<\dx, dy\ — DjD- {x, y}, the module of the delta functions 
with the support (0,0). Then we compute 

15 : DExt(M,N) 

05 = HashTable{ => QQ° } 

1 => qq' 

2 => qQ= 
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As before, once we know the dimension of Hom£)(Af, iV), we can compute 
the solutions of M in A'' by a brute force method. 

1. For given holonomic systems AI and N, evaluate the dimension d of 
HomD(Af, N) by the homological duahty method. 

2. Filter N by finite-dimensional vector spaces F^{N) and search for solutions 
in Fi (N) for increasing i until d linearly independent solutions are found. 

For instance in step 2, it N ~ D/J, then we can use the induced Bernstein 
filtration B where B^{D / J) consists of residues of elements L £ D whose total 
degree is less than or equal to i. 
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